#include <iostream>
#include <conio.h>
#include <cmath>

using namespace std;

double const M = 5.97219e24;	//units: kg
double const G = 6.67e-11;		//units: Nm2/kg2

void calculate_new_values(double* R, double* V, double* A, double m, double F, double dt, double* new_R, double* new_V, double* new_A);

int main(){
	double dt, m, F;
	

	double r[3] = {4082000.0, 4082000.0, 4082000.0}, new_r[3]= {0.0, 0.0, 0.0};		//Initial position at a 700 km orbit
	double v[3] = {0.0, 10000.0, 0.0}, new_v[3]= {0.0, 0.0, 0.0};					//Initial velocity 10 km/s
	double a[3] = {0.0, 0.0, 0.0}, new_a[3]={0.0, 0.0, 0.0};


	cout << "Initial Position:\n" << "\nX:" <<r[0] << "\nY:" <<r[1] <<"\nZ:" <<r[2] << endl << endl;
	cout << "Initial Velocity:\n" << "\nX:" <<v[0] << "\nY:" <<v[1] <<"\nZ:" <<v[2] << endl << endl;
	cout << "Initial Acceleration:\n" << "\nX:" <<a[0] << "\nY:" <<a[1] <<"\nZ:" <<a[2] << endl;

	cout << "\nDefine the time step in seconds: ";
	cin >> dt;

	cout << "\nDefine the mass of the spacecraft in kilogram: ";
	cin >> m;

	cout << "\nDefine the thrust in Newton: ";
	cin >> F;

	for (int i=0; i<10; i++){
		
		calculate_new_values(r,v,a,m,F,dt,new_r,new_v,new_a);

		cout << "New Position:\n" << "\nX:" <<new_r[0] << "\nY:" <<new_r[1] <<"\nZ:" <<new_r[2] << endl << endl;
		cout << "New Velocity:\n" << "\nX:" <<new_v[0] << "\nY:" <<new_v[1] <<"\nZ:" <<new_v[2] << endl << endl;
		cout << "New Acceleration:\n" << "\nX:" <<new_a[0] << "\nY:" <<new_a[1] <<"\nZ:" <<new_a[2] << endl << endl;
		
		for (int j=0; j<3; j++){
			r[j]=new_r[j];
			v[j]=new_v[j];
			a[j]=new_a[j];
		}
	}




	_getch();
	return 0;

}




void calculate_new_values(double* R, double* V, double* A, double m, double F, double dt, double* new_R, double* new_V, double* new_A){

	
	//R = vector that holds state variables
	//dt = timestep in seconds
	//N = number of components in vector R
	//Rd = vector that holds the derivatives of the state variables

	R[0];			//Position X
	R[1];			//Position Y
	R[2];			//Position Z
	V[0];			//Velocity X
	V[1];			//Velocity Y
	V[2];			//Velocity Z
	

	double r2 = R[0]*R[0] + R[1]*R[1] + R[2]*R[2];		//r^2
	double r = sqrt(r2);
	double r3 = r2*r;

	double Formula = ((G*M*m)+(F*r2))/(m*r3);

	for (int i = 0; i<3; i++){
		new_A[i] = (-R[i] * Formula);
		new_V[i] = (new_A[i] * dt) + V[i];
		new_R[i] = (new_V[i] * dt) + R[i];
	}

 

}